---
title: "Hide and Seek: The connection between false beliefs and perceptions of government transparency"
date: "2022-03-13"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

## Set the working directory and load the packages

```{r}
# Set working directory
setwd(dirname(rstudioapi::getActiveDocumentContext()$path))

# Load packages
#library(pacman)
pacman::p_load(tidyverse, haven, psych, scales, sandwich, gridExtra, irr, nFactors, kableExtra, estimatr, jtools, ggeffects, ggstance, nnet, MNLpred, texreg, stargazer, ebal, survey, weights, survey, broom.mixed, patchwork, WeightIt, cobalt, poLCA, GPArotation, apaTables, lemon, specr)
```

## Create folders to save the figures and tables

```{r}
if (!file.exists("Figures")) dir.create("Figures")
if (!file.exists("Tables")) dir.create("Tables")
```

## Functions

```{r}
# To rescale variables on 0-1 scales
range01 <- function(x, ...){(x - min(x, ...)) / (max(x, ...) - min(x, ...))}

# To avoid repeating all these arguments in stargazer when creating the tables
stargazer2 <- purrr::partial(stargazer, type="html", model.names = F,
          title = "", label = "", dep.var.caption = "", dep.var.labels.include = FALSE, omit.stat = "f",
          digits = 2, df = FALSE, omit=to_omit, star.cutoffs = c(0.1, 0.05, 0.01, 0.001),
          star.char = c("+", "*", "**", "***"), 
          notes = c("+ p<0.1; * p<0.05; ** p<0.01; *** p<0.001"), 
          notes.append = F)
to_omit <- c("Constant", "regionPrairies", "regionBC", "regionMaritimes", "regionQuebec")

# To create a dataset of weighted means, se, and variable names
wt_prop_se <- function(x, design){
  require(tidyverse)
  require(survey)
  x %>% 
  gather(key = "key", value = "value") %>% 
  dplyr::group_by(key) %>% 
  dplyr::summarize(mean = as.data.frame(svymean(~value, design, na.rm = T))[,1],
                   SE = as.data.frame(svymean(~value, design, na.rm = T))[,2],
                   cat = str_remove(rownames(as.data.frame(svymean(~value, design, na.rm = T))), "value"))
}
```

## Load data

```{r}
dat <- read_csv("./HKS_covid_canada.csv")
```

## Recode variables for main analyses

#### Dependent variables (perceptions of transparency)

```{r}
### Lacks transparency
dat <- dat %>% 
  mutate(transparent = dplyr::recode(as.numeric(Q74r1_V3), `1`=1, `2`=0.67, `3`=0.33, `4`=0, `9`=NA_real_))

### Hides information
dat <- dat %>% 
  mutate(gov_hiding10_V3 = ifelse(Q36_V3 != 99, Q36_V3, NA),
         gov_hiding_V3 = gov_hiding10_V3/10,
         gov_hiding10_V2 = ifelse(Q36_V2 != 99, Q36_V2, NA),
         gov_hiding_V2 = gov_hiding10_V2/10,
         gov_hiding_var23 =  gov_hiding10_V3 - gov_hiding10_V2,
         gov_hiding10_V1 = ifelse(Q36_V1 != 99, Q36_V1, NA),
         gov_hiding_V1 = gov_hiding10_V1/10)
```

#### Other dependent variables (satisfaction with democracy, vote intentions)

```{r}
dat <- dat %>% 
  ### Satisfaction with democracy
  mutate(satisf_democ_V3 = dplyr::recode(as.numeric(Q46_V3),
                                         `1`=1, `2`=0.67, `3`=0.33, `4`=0, `9`=NA_real_),
         satisf_democ_V2 = dplyr::recode(as.numeric(Q46_V2),
                                         `1`=1, `2`=0.67, `3`=0.33, `4`=0, `9`=NA_real_),
         satisf_democ_var23 = satisf_democ_V3 - satisf_democ_V2,
  ### Vote intentions for the governing federal party
         vote_LPC_V3 = ifelse(Q30_V3 == 1, 1, 0),
         vote_LPC_V2 = case_when(
           Q30_V2 == 1 ~ 1,
           Q30_V2 %in% 2:99 ~ 0,
           TRUE ~ NA_real_),
         vote_LPC_var23 = vote_LPC_V3 - vote_LPC_V2,
  ### Vote intentions for the governing provincial party
         vote_prov_V3 = case_when(
           Q30AB_V3 == 4 | Q30BC_V3 == 3 | Q30PEI_V3 == 4 |
             Q30MB_V3 == 3 | Q30NB_V3 == 2 | Q30NS_V3 == 1 |
             Q30ON_V3 == 1 | Q30QC_V3 == 1 | Q30SK_V3 == 5 | 
             Q30NF_V3 == 1 ~ 1,
           Q30AB_V3 %in% c(1,2,3,5:99) | Q30BC_V3 %in% c(1,2, 4:99) | 
             Q30PEI_V3 %in% c(1,2,3, 5:99) |
             Q30MB_V3 %in% c(1,2, 4:99) | Q30NB_V3 %in% c(1, 3:99) | Q30NS_V3 %in% c(2:99) |
             Q30ON_V3 %in% c(2:99) | Q30QC_V3 %in% c(2:99) | Q30SK_V3 %in% c(1,2,3,4, 6:99) | 
             Q30NF_V3 %in% c(2:99) ~ 0,
           TRUE ~ NA_real_),
         vote_prov_V2 = case_when(
           Q30AB_V2 == 4 | Q30BC_V2 == 3 | Q30PEI_V2 == 4 |
             Q30MB_V2 == 3 | Q30NB_V2 == 2 | Q30NS_V2 == 1 |
             Q30ON_V2 == 1 | Q30QC_V2 == 1 | Q30SK_V2 == 5 | 
             Q30NF_V2 == 1 ~ 1,
           Q30AB_V2 %in% c(1,2,3,5:99) | Q30BC_V2 %in% c(1,2, 4:99) | 
             Q30PEI_V2 %in% c(1,2,3, 5:99) |
             Q30MB_V2 %in% c(1,2, 4:99) | Q30NB_V2 %in% c(1, 3:99) | Q30NS_V2 %in% c(2:99) |
             Q30ON_V2 %in% c(2:99) | Q30QC_V2 %in% c(2:99) | Q30SK_V2 %in% c(1,2,3,4, 6:99) | 
             Q30NF_V2 %in% c(2:99) ~ 0,
           TRUE ~ NA_real_),
         vote_prov_var23 = vote_prov_V3 - vote_prov_V2)
```

#### False beliefs

```{r}
dat <- dat %>%
  ### Average truthfulness of the false statements
  mutate(misinfo = dat$MISINFOr1_V3 + dat$MISINFOr2_V3 + 
           dat$MISINFOr3_V3 + dat$MISINFOr6_V3 + dat$MISINFOr7_V3 +
           dat$MISINFOr9_V3,
         misinfo10 = misinfo/6,
         misinfo = misinfo/6/10,
  ### Average truthfulness of the true statements
         true_all = MISINFOr4_V3 + MISINFOr5_V3 + MISINFOr8_V3,
         true_all = true_all/3/10,
  ### Difference in perceived truthfulness of false and true statements
         susceptibility = misinfo - true_all,
  ### Rename false beliefs
        restrict_rights = MISINFOr1_V3,
        power_china = MISINFOr2_V3,
        benefit_companies = MISINFOr3_V3,
        co2_masks = MISINFOr6_V3,
        hydroxychloroquine = MISINFOr7_V3,
        false_positives = dat$MISINFOr9_V3,
  ### Number of false beliefs
         restrict_rights01 = ifelse(dat$MISINFOr1_V3 > 5, 1, 0),
         power_china01 = ifelse(MISINFOr2_V3 > 5, 1, 0),
         benefit_companies01 = ifelse(MISINFOr3_V3 > 5, 1, 0),
         co2_masks01 = ifelse(MISINFOr6_V3 > 5, 1, 0),
         hydroxychloroquine01 = ifelse(MISINFOr7_V3 > 5, 1, 0),
         false_positives01 = ifelse(MISINFOr9_V3 > 5, 1, 0),
         nmisinfo = restrict_rights01 + power_china01 + benefit_companies01 + 
           co2_masks01 + hydroxychloroquine01 + false_positives01,
  ### Number of false beliefs as a categorical variable (for Figure 2)
         nmisinfo_cat = case_when( 
           nmisinfo == 0 ~ "None (57%)",
           nmisinfo %in% c(1,2) ~ "One or two (20%)",
           nmisinfo %in% c(3,4) ~ "Three or four (12%)",
           nmisinfo %in% c(5,6) ~ "Five or six (10%)",
           TRUE ~ NA_character_),
  ### True vs. false vs. unsure (for Figure A2) 
        restrict_rights2 = ifelse(MISINFOr1_V3 == 5, 0.5, restrict_rights01),
        power_china2 = ifelse(MISINFOr2_V3 == 5, 0.5, power_china01),
        benefit_companies2 = ifelse(MISINFOr3_V3 == 5, 0.5, benefit_companies01),
        co2_masks2 = ifelse(MISINFOr6_V3 == 5, 0.5, co2_masks01),
        hydroxychloroquine2 = ifelse(MISINFOr7_V3 == 5, 0.5, hydroxychloroquine01),
        false_positives2 = ifelse(MISINFOr9_V3 == 5, 0.5, false_positives01),
  ### Holding false beliefs (dummy)
  #### For showing that the results hold when using ebal in Appendix C
        misinfo01 = ifelse(nmisinfo >= 1, 1, 0),
  ### Individual false beliefs on 0-1 scales (Table C2 and C3)
        restrict_rights_scaled = restrict_rights/10,
        power_china_scaled = power_china/10,
        benefit_companies_scaled = benefit_companies/10,
        co2_masks_scaled = co2_masks/10,
        hydroxychloroquine_scaled = hydroxychloroquine/10,
        false_positives_scaled = false_positives/10)
```

#### Control variables

```{r}
# Trust
dat <- dat %>% 
  ### Trust in the federal government
  mutate(trust_fed_V3 = dplyr::recode(as.numeric(Q7r4_V3), 
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_),
         trust_fed_V2 = dplyr::recode(as.numeric(Q7r4_V2),
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_),
         trust_fed_var23 = trust_fed_V3 - trust_fed_V2,
  ### Trust in the provincial government
         trust_prov_V3 = dplyr::recode(as.numeric(Q7r5_V3),
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_),
         trust_prov_V2 = dplyr::recode(as.numeric(Q7r5_V2), 
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_),
         trust_prov_var23 = trust_prov_V3 - trust_prov_V2,
  ### Trust in scientists
         trust_scien_V3 = dplyr::recode(as.numeric(Q7r6_V3),
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_), 
         trust_scien_V2 = dplyr::recode(as.numeric(Q7r6_V2), 
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_), 
  ### Trust in journalists
         trust_journ_V3 = dplyr::recode(as.numeric(Q7r9_V3),
                          `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, `9` = NA_real_), 
  ### Generalized social trust
         general_trust_V3 = ifelse(Q4_V3 == 1, 1, 0))

# Partisan identification
### Federal PID
dat <- dat %>% 
  ### LPC is the governing federal party
  mutate(LPC_V3 = case_when( 
  Q42_V3 == 1 ~ 1,
  Q42_V3 %in% c(2:99) ~ 0,
  Q41B_V3 %in% c(2, 9) ~ 0,
  Q41_V3 == 9 ~ 0,
  TRUE ~ NA_real_),
  ### Identification with the governing party at t-1
  LPC_V2 = case_when(
  Q42_V2 == 1 ~ 1,
  Q42_V2 %in% c(2:99) ~ 0,
  Q41_V2 %in% c(2, 9) ~ 0,
  TRUE ~ NA_real_),
  ### CPC stands for Conservative Party of Canada
  CPC_V3 = case_when(
  Q42_V3 == 2 ~ 1,
  Q42_V3 %in% c(1, 3:99) ~ 0,
  Q41B_V3 %in% c(2, 9) ~ 0,
  Q41_V3 == 9 ~ 0,
  TRUE ~ NA_real_),
  ### NDP stands for New Democratic Party
  NDP_V3 = case_when(
  Q42_V3 == 3 ~ 1,
  Q42_V3 %in% c(1:2, 4:99) ~ 0,
  Q41B_V3 %in% c(2, 9) ~ 0,
  Q41_V3 == 9 ~ 0,
  TRUE ~ NA_real_))

### Provincial PID (dummy for identification with the governing party)
dat <- dat %>% 
  mutate(PID_prov_V3  = case_when(
    Q44AB_V3 == 4 | Q44BC_V3 == 3 | Q44PEI_V3 == 4 |
      Q44MB_V3 == 3 | Q44NB_V3 == 2 | Q44NS_V3 == 1 |
      Q44ON_V3 == 1 | Q44QC_V3 == 1 | Q44SK_V3 == 5 | 
      Q44NF_V3 == 1 ~ 1,
    TRUE ~ 0),
  PID_prov_V2 = case_when(
    Q44AB_V2 == 4 | Q44BC_V2 == 3 | Q44PEI_V2 == 4 |
      Q44MB_V2 == 3 | Q44NB_V2 == 2 | Q44NS_V2 == 1 |
      Q44ON_V2 == 1 | Q44QC_V2 == 1 | Q44SK_V2 == 5 | 
      Q44NF_V2 == 1 ~ 1,
    Q43_V2 %in% c(2,9) ~ 0,
    !is.na(Q44AB_V2) | !is.na(Q44BC_V2) | !is.na(Q44PEI_V2) |
      !is.na(Q44MB_V2) | !is.na(Q44NB_V2) | !is.na(Q44NS_V2) |
      !is.na(Q44ON_V2) | !is.na(Q44QC_V2) | !is.na(Q44SK_V2) | 
      !is.na(Q44NF_V2) ~ 0,
    TRUE ~ NA_real_),
  PID_prov_var23 = PID_prov_V3 - PID_prov_V2)

# Ideology (rescaled from 0 to 1)
dat$ideol_V3 <- ifelse(dat$Q40_V3 != 99, dat$Q40_V3/10, NA)
dat$ideol_V2 <- ifelse(dat$Q40_V2 != 99, dat$Q40_V2/10, NA)

# Satisfaction with the handling of the pandemic
dat <- dat %>% 
  mutate(handling_fed_V3 = dplyr::recode(as.numeric(Q14A_V3), 
                  `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, .default = NA_real_),
         handling_fed_V2 = dplyr::recode(as.numeric(Q14A_V2), 
                  `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, .default = NA_real_),
         handling_prov_V3 = dplyr::recode(as.numeric(Q14B_V3), 
                  `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, .default = NA_real_),
         handling_prov_V2 = dplyr::recode(as.numeric(Q14B_V2), 
                  `1` = 1, `2` = 0.67, `3` = 0.33, `4` = 0, .default = NA_real_),
         handling_fed_var23 = handling_fed_V3 - handling_fed_V2,
         handling_prov_var23 = handling_prov_V3 - handling_prov_V2)

# Media consumption 
dat$soc_media <- dat$COROINFOr4_V3/5
dat$tv <- dat$COROINFOr1_V3/5
dat$prints <- dat$COROINFOr3_V3/5
dat$radio <- dat$COROINFOr2_V3/5

# Socio-demographics
dat <- dat %>% 
  mutate(age = range01(dat$age_V3, na.rm = T),
         female = range01(dat$sexe_V3, na.rm = T),
         educ = case_when(
           SCOLA_V3 %in% c(1, 2) ~ 0,
           SCOLA_V3 == 3 ~ 0.5,
           SCOLA_V3 == 4 ~ 1,
           TRUE ~ NA_real_))

# Region
dat <- dat %>% 
  mutate(region = case_when(
    PROV2_V3 == 1 ~ "BC", 
    PROV2_V3 %in% c(2,3) ~ "Prairies",
    PROV2_V3 == 4 ~ "Ontario",
    PROV2_V3 == 5 ~ "Quebec",
    PROV2_V3 == 6 ~ "Maritimes",
    TRUE ~ NA_character_),
    region = relevel(factor(region), ref = "Ontario"),
### Regions as dummies (for entropy balancing and predictions)
    Ontario = ifelse(region == "Ontario", 1, 0),
    Quebec = ifelse(region == "Quebec", 1, 0),
    Maritimes = ifelse(region == "Maritimes", 1, 0),
    BC = ifelse(region == "BC", 1, 0),
    Prairies = ifelse(region == "Prairies", 1, 0))

# Life satisfaction
dat$satisf_life <- ifelse(dat$Q1_V3 != 99, dat$Q1_V3/10, NA)

# Emotions 
dat$fear <- ifelse(dat$Q11r1_V3 != 99, dat$Q11r1_V3/10, NA)
dat$anger <- ifelse(dat$Q11r3_V3 != 99, dat$Q11r3_V3/10, NA)

# Save design
dat_design <- svydesign(ids=~1, data=dat, weights = dat$pond_V3)
```

## Main models

#### Figure 1: False beliefs and perceptions of transparency

```{r}
controls <- "trust_scien_V3 + LPC_V3 + PID_prov_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + region"

# Estimate linear regression

## Perceptions that governments hide information
### Standard SEs
hiding_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ misinfo +", controls)))
### With robust standard errors
hiding_lm_robust <- lm_robust(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ misinfo +", controls)))
summary(hiding_lm_robust)

## Perceptions that governments hide information
### Standard SEs
transparency_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("transparent ~ misinfo +", controls)))
### With robust standard errors
transparency_lm_robust <- lm_robust(data = dat, weights = pond_V3, 
   formula = formula(paste0("transparent ~ misinfo +", controls)))
summary(transparency_lm_robust)

# Plot unstandardized coefficients
figure_1 <- plot_summs(hiding_lm, transparency_lm, 
           inner_ci_level = .9,
           scale = F, robust = T,
           model.names = c( "Hides information", "Lacks transparency"),
           colors = "Qual1",
           coefs = c("False beliefs (0-10)" = "misinfo",
                     "Trust scientists (1-4)" = "trust_scien_V3",
                     "Federal PID (0-1)" = "LPC_V3",
                     "Provincial PID (0-1)" = "PID_prov_V3",
                     "Left-right ideology (0-10)" = "ideol_V3",
                     "Social media news consumption (1-5)" = "soc_media",
                     "TV news consumption (1-5)" = "tv",
                     "Print news consumption (1-5)" = "prints",
                     "Age (1-7)" = "age",
                     "Female (0-1)" = "female",
                     "Education (1-4)" = "educ"))+
  scale_x_continuous(limits = c(-0.18, 0.80), 
                     breaks = c(-0.2, 0, 0.2, 0.4, 0.6, 0.8))+
  theme_minimal()+
  labs(y = "",
       x = "Effect on outcome (0-1)")+
  guides(shape = "none")+
  theme(axis.text=element_text(size=10),
        axis.title=element_text(size=12),
        legend.text = element_text(size=10)); figure_1

ggsave(figure_1, file = "./Figures/figure_1.png", height = 5, width = 8, dpi = 500)
```

#### Figure 2: Average perceived likelihood (0-10) that governments are hiding information about the pandemic in each survey wave based on the number of false statements endorsed

```{r}
# Estimate weighted means for each survey wave
false_beliefs_W1 <- svyby(~gov_hiding10_V1, ~nmisinfo_cat, dat_design, svymean, na.rm = T)
false_beliefs_W2 <- svyby(~gov_hiding10_V2, ~nmisinfo_cat, dat_design, svymean, na.rm = T)
false_beliefs_W3 <- svyby(~gov_hiding10_V3, ~nmisinfo_cat, dat_design, svymean, na.rm = T)

# Prepare data for graph
false_beliefs_W1 <- false_beliefs_W1 %>% 
  mutate(Wave = "1") %>% 
  rename("gov_hiding" = "gov_hiding10_V1")
false_beliefs_W2 <- false_beliefs_W2 %>% 
  mutate(Wave = "2") %>% 
  rename("gov_hiding" = "gov_hiding10_V2")
false_beliefs_W3 <- false_beliefs_W3 %>% 
  mutate(Wave = "3") %>% 
  rename("gov_hiding" = "gov_hiding10_V3")
false_beliefs_W <- rbind(false_beliefs_W1, false_beliefs_W2, false_beliefs_W3)
false_beliefs_W <- false_beliefs_W %>% 
  mutate(nmisinfo_cat = factor(nmisinfo_cat, levels = c("None (57%)", "One or two (20%)", "Three or four (12%)", "Five or six (10%)")))

# Make graph
hides_by_beliefs_plot <- ggplot(false_beliefs_W, aes(x = Wave, y = gov_hiding, 
                            ymin = gov_hiding-1.96*se, ymax = gov_hiding+1.96*se, 
                            col = nmisinfo_cat))+
  geom_point(size = 2, position = position_dodge(0.2))+
  geom_errorbar(width = 0.01, position = position_dodge(0.2))+
  scale_color_manual(values=rev(RColorBrewer::brewer.pal(4, "Set1")))+
  labs(y = "Hides information (0-10)",
       col = "Number of false\nstatements endorsed")+
  theme_minimal()+
   theme(axis.text=element_text(size=11),
        axis.title=element_text(size=12)); hides_by_beliefs_plot

ggsave(hides_by_beliefs_plot, file = "./Figures/figure_2.png", height = 5, dpi = 500)
```

#### Linear models: False beliefs and changes in perceptions of transparency

```{r}
# Three model specifications
controls_var <- "trust_scien_V2 + LPC_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + gov_hiding10_V2 + region"

controls_var1 <- "handling_fed_var23 + handling_prov_var23 + trust_scien_V2 + LPC_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + gov_hiding10_V2 + region"

controls_var2 <- "trust_fed_var23 + trust_prov_var23 + handling_fed_var23 + handling_prov_var23 + trust_scien_V2 + LPC_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + gov_hiding10_V2 + region"

# Estimate linear regression models
hiding_var_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ misinfo +", controls_var)))
hiding_var_lm_robust <- lm_robust(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ misinfo +", controls_var)))
summary(hiding_var_lm_robust)

## Estimate the model with the false beliefs variable to compare the model fit
hiding_var_lm_nomisinfo <- lm_robust(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~", controls_var)))
summary(hiding_var_lm_nomisinfo)

# Alternative specifications (see Table B4)
hiding_var_lm1 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ misinfo +", controls_var1)))

hiding_var_lm2 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ misinfo +", controls_var2)))
```

### Figure 3: False beliefs and justification for transparency perceptions (open-ended)

```{r}
# Calculate percentage of agreement
dat <- dat %>% 
  mutate(intercoder = ifelse(open_ended1 == open_ended2, 1, 0))
answered <- dat %>% 
  filter(Q36Br1_V3 != "" & !is.na(Q36Br1_V3))
mean(answered$intercoder, na.rm = T)

# Calculate Cohen's Kappa
kappa2(answered[, c("open_ended1", "open_ended2")], "squared")

# Deal with observations for which the coders do not agree
## Remove these observations
dat$open_ended <- ifelse(dat$open_ended1 == dat$open_ended2, dat$open_ended1, NA)
# Randomly select between coder 1 or 2
set.seed(123)
dat$rnum <- sample(c(1,2), length(dat$open_ended), replace = T)
dat$open_ended5 <- ifelse(dat$rnum == 1, dat$open_ended1, dat$open_ended2)

# Keep in categories of interest
dat_open <- dat %>% 
  mutate(open_ended = ifelse(open_ended %in% c("positive", "handling", "mistrust", "manipulation"), open_ended, NA),
         open_ended1 = ifelse(open_ended1 %in% c("positive", "handling", "mistrust", "manipulation"), open_ended1, NA),
         open_ended2 = ifelse(open_ended2 %in% c("positive", "handling", "mistrust", "manipulation"), open_ended2, NA),
         open_ended5 = ifelse(open_ended5 %in% c("positive", "handling", "mistrust", "manipulation"), open_ended5, NA),
         open_ended = factor(open_ended, levels = c("positive", "handling", "mistrust", "manipulation")),
         open_ended1 = factor(open_ended1, levels = c("positive", "handling", "mistrust", "manipulation")),
         open_ended2 = factor(open_ended2, levels = c("positive", "handling", "mistrust", "manipulation")),
         open_ended5 = factor(open_ended5, levels = c("positive", "handling", "mistrust", "manipulation")))

# Estimate multinomial model 
misinfo_multinom5 <- multinom(open_ended5 ~ misinfo + trust_scien_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + Maritimes + Quebec + Prairies + BC, data = dat_open, maxit = 200, weights = pond_V3, Hess = T)

# Generate predictions
pred5 <- mnl_pred_ova(model = misinfo_multinom5,
                      data = dat_open,
                      x = "misinfo",
                      by = 0.01666667,
                      seed = 123, 
                      nsim = 1000, 
                      probs = c(0.025, 0.975)) 

# Calculate first difference in predicted probability (as reported in the text)
mnl_fd2_ova(model = misinfo_multinom5,
                     data = dat_open,
                     x = "misinfo",
                     seed = 123,
                     value1 = min(dat_open$misinfo),
                     value2 = max(dat_open$misinfo),
                     nsim = 1000)

# Create data for the plot
pred5$plotdata <- pred5$plotdata %>% 
  mutate(level = factor(open_ended5, levels = c("manipulation", "mistrust", "handling", "positive"),
                        labels = c("Manipulation", "Mistrust", "Handling", "Positive")))

# Plot the results
figure_3 <- ggplot(data = pred5$plotdata, aes(x = misinfo, 
                                  y = mean,
                                  ymin = lower, ymax = upper)) +
  geom_ribbon(aes(fill = level, col = level), alpha = 0.1) + 
  geom_line(aes(col = level)) +
  facet_wrap(.~ level, ncol = 2) +
  scale_color_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_fill_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  scale_x_continuous(breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  theme_bw() +
  labs(y = "Probability of falling into each category",
       x = "Average perceived truthfulness of the false statements")+
  guides(fill = "none", col = "none")+
  theme_minimal()+
  theme(axis.text=element_text(size=12),
        axis.title=element_text(size=13),
        legend.text = element_text(size = 11),
        strip.text = element_text(size = 13));figure_3

ggsave(figure_3, file = "./Figures/figure_3.png", height = 5, width = 8, dpi = 500)
```

## Appendices

## Appendix A: Descriptive statistics

#### Table A1: Descriptive statistics table

```{r}
select_desc <- dplyr::select(dat, transparent, gov_hiding_V3, gov_hiding_V2, gov_hiding_var23, misinfo, vote_LPC_V3, vote_LPC_V2, vote_LPC_var23, vote_prov_V3, vote_prov_V2, vote_prov_var23, trust_fed_V3, trust_fed_V2, trust_fed_var23, trust_prov_V3, trust_prov_V2, trust_prov_var23, satisf_democ_V3, satisf_democ_V2, satisf_democ_var23, handling_fed_V3, handling_fed_V2, handling_fed_var23, handling_prov_V3, handling_prov_V2, handling_prov_var23, trust_scien_V3, trust_scien_V2, LPC_V3, PID_prov_V3, ideol_V3, soc_media, tv, prints, age, female, educ)

varnames_desc <- c("Lacks transparency", "Hides information", "Hides information (t-1)", "Δ Hides information", "False beliefs",  "Vote federal", "Vote federal (t-1)", "Δ Vote federal", "Vote provincial", "Vote provincial (t-1)", "Δ Vote provincial", "Trust federal",  "Trust federal (t-1)", "Δ Trust federal", "Trust provincial", "Trust provincial (t-1)",  "Δ Trust provincial", "Satisfaction democracy", "Satisfaction democracy (t-1)", "Δ Satisfaction democracy", "Handling federal", "Handling federal (t-1)", "Δ Handling federal", "Handling provincial", "Handling provincial (t-1)", "Δ Handling provincial", "Trust scientists", "Trust scientists (t-1)", "Federal PID", "Provincial PID", "Ideology", "Social media news consumption", "TV news consumption", "Print news consumption", "Age", "Female", "Education")

desc_transparency <- as.data.frame(psych::describe(select_desc, fast = T))
rownames(desc_transparency) <- varnames_desc
desc_transparency <- desc_transparency %>% 
  round(digits = 3) %>% 
  dplyr::select(-vars, -range, -se)  

desc_transparency %>%  
kbl() %>% 
  kable_classic(html_font = "Calibri") %>% 
  save_kable(file = "./Tables/table_A1.html", self_contained = T)
```

#### Figure A1: Distribution of the perceived truthfulness of each false statement

```{r}
figure_A1 <- dat %>%
  dplyr::select(restrict_rights, power_china, benefit_companies, 
         co2_masks, hydroxychloroquine, false_positives) %>% 
  sapply(., as.factor) %>% 
  as.data.frame() %>% 
  wt_prop_se(., dat_design) %>% 
  mutate(key = factor(key, levels = rev(c("false_positives", "hydroxychloroquine",
                                        "co2_masks", "benefit_companies", "power_china",
                                        "restrict_rights")), 
                       labels = rev(c("False\npositives", "Hydroxy\nchloroquine", "CO2\nMasks", 
                                  "Benefit\ncompanies", "Power\nChina", "Restrict\nrights"))),
         cat = factor(cat, levels = 0:10)) %>% 
  ggplot(., aes(x = cat, y=mean, 
                ymin=mean-1.96*SE, ymax=mean+1.96*SE, fill = key))+
  geom_bar(stat = "identity")+
  geom_errorbar(aes(col = key), width=0.1, alpha=0.9, size=0.3, alpha = 0.5)+
  facet_wrap(~key, ncol = 1, strip.position = "right")+
  scale_x_discrete(breaks = c("0","1","2","3","4","5","6","7","8","9","10"), labels = c("0\n", "1\nNot likely at all to be true", "2\n", "3\n", "4\n", "5\n", "6\n", "7\n", "8\n", "9\nExtremely likely to be true", "10\n"))+
  scale_y_continuous(labels = percent_format(accuracy = 1))+
  scale_fill_manual(values = RColorBrewer::brewer.pal(5, "Set1")[1:6])+
  scale_color_manual(values = c("red", "cyan4", "green4", "purple", "orange", "grey"))+
  guides(fill = "none", col = "none")+
  labs(x = "Perceived truthfulness (0-10)", 
       y = "Percentage of respondents")+
  theme_minimal()+
  theme(panel.spacing = unit(0.25, "lines"));figure_A1

ggsave(figure_A1, file = "./Figures/figure_A1.png", width = 5, height = 5)
```

#### Figure A2: Weighted percentage of respondents endorsing each false statement

```{r}
figure_A2 <- dat %>% 
  dplyr::select(restrict_rights2, power_china2, benefit_companies2, co2_masks2, hydroxychloroquine2, false_positives2, pond_V3) %>% 
  pivot_longer(cols = c("restrict_rights2", "power_china2", "benefit_companies2", "co2_masks2", "hydroxychloroquine2", "false_positives2"),names_to = "names", values_to = "value") %>% 
  mutate(names = rep(c("Restrict rights", "Power China", "Benefit companies", "CO2 masks", "Hydroxychloroquine", "False positives"), 1005)) %>% 
 # mutate(names2 = factor(names, labels = rev(c("Restrict rights", "Power China", "Benefit companies", "CO2 masks", "Hydroxychloroquine", "False positives")))) %>% 
  ggplot(aes(x = reorder(names, -value), fill = factor(value)))+
  geom_bar(aes(y=..count../sum(..count..)*100, weight = pond_V3), position = "fill")+
  labs(x = "",
       y = "Percentage of respondents",
       fill = "Truthfulness")+
  scale_fill_manual(labels = c("False", "Unsure", "True"), values = RColorBrewer::brewer.pal(4, "Set1")[3:1])+
   scale_y_continuous(breaks = c(0, 0.20, 0.40, 0.60, 0.80, 1), 
                      labels = percent_format(accuracy = 1))+
  theme_minimal()+
  theme(axis.text=element_text(size=12),
        axis.title=element_text(size=13),
        legend.text=element_text(size=11))+
  coord_flip();figure_A2

ggsave(figure_A2, file = "./Figures/figure_A2.png", height = 5)
```

### Figure A3: Number of false statements endorsed

```{r}
figure_A3 <- as.data.frame(svymean(~factor(nmisinfo), dat_design, na.rm = T)) %>% 
  ggplot(aes(x = 0:6, y = mean, ymin = mean-1.96*SE, ymax = mean+1.96*SE)) +
  geom_bar(stat = "identity", fill = "darkblue")+
  geom_errorbar(width=0.1, alpha=0.9, size=0.7, col = "blue2")+
  labs(x = "Number of false statements endorsed", 
       y = "Percentage of respondents")+
  scale_x_continuous(breaks = c(0:6))+
  scale_y_continuous(breaks = c(0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6), 
                     labels = percent_format(accuracy = 1))+
  theme_minimal()+
  theme(axis.text=element_text(size=12),
        axis.title=element_text(size=13));figure_A3

ggsave(figure_A3, file = "./Figures/figure_A3.png", height = 5)
```

### Table A2: Correlation between the perceived truthfulness of false statements

```{r}
# Select variables
misinfo_select <- dplyr::select(dat, restrict_rights, power_china, benefit_companies, co2_masks, hydroxychloroquine, false_positives) 
# Rename the columns
colnames(misinfo_select) <- c("Restrict rights", "Power China", "Benefit companies", "CO2 intoxication", "Hydroxychloroquine", "False positives")

# Create the correlation matrix
mcor<-round(cor(misinfo_select), 2)
upper<-mcor
upper[upper.tri(mcor, diag = TRUE)]<-""
upper <- upper %>% 
  as.data.frame()

# Create the table
upper %>% 
  kbl() %>% 
  kable_classic(html_font = "Calibri") %>% 
  save_kable(file = "./Tables/table_A2.html", self_contained = T)

# Correlation tests (reported in the main text)
cor.test(dat$hydroxychloroquine, dat$power_china, use = "complete.obs")
cor.test(dat$false_positives, dat$restrict_rights, use = "complete.obs")
```

## Other statistics discussed in the text 

#### Principal component analysis

```{r}
# Select false beliefs variables
misinfo_df <- dplyr::select(dat, restrict_rights, power_china, benefit_companies, co2_masks, hydroxychloroquine, false_positives)

# Principal component analysis
fit <- princomp(misinfo_df, cor=T)
vaf <- summary(fit) # print variance accounted for
pc_loadings <- loadings(fit) # pc loadings
plot(fit,type="lines") # scree plot
```

#### Cronbach alpha

```{r}
# Calculate cronbach alpha
psych::alpha(misinfo_df, check.keys = T)
# Cronbach alpha is 0.91
```

#### Figure A4: Distribution of perceptions of transparency 

```{r}
# Lacks transparency (Left pane)
mean_transparent <- svymean(~factor(transparent), dat_design, na.rm = T)
confint_transparent <- confint(svymean(~factor(transparent), dat_design, na.rm = T))

plot_transparent <- data.frame(lack_transparency = factor(c("Completely disagree", "Somewhat disagree", 
                              "Somewhat agree", "Completely agree"), 
                              levels = c("Completely disagree", "Somewhat disagree", 
                              "Somewhat agree", "Completely agree")),
                              prop = c(mean_transparent[[1]], mean_transparent[[2]],
                                       mean_transparent[[3]], mean_transparent[[4]]),
                              lower = confint_transparent[1:4],
                              upper = confint_transparent[5:8]) %>% 
  ggplot(., aes(x = lack_transparency, y = prop, ymin = lower, ymax = upper)) +
  geom_bar(stat = "identity", fill = "darkblue")+
  geom_errorbar(width=0.1, alpha=0.9, size=0.7, col = "blue2")+
  labs(x = "\n Provincial government lacks transparency \n on what drives its decisions",
       y = "Percentage of respondents")+
  scale_y_continuous(labels = percent_format(accuracy = 1))+
  scale_x_discrete(labels = c("Completely \n disagree", "Somewhat \n disagree", 
                              "Somewhat \n agree", "Completely \n agree"))+
  theme_minimal()+
  theme(axis.text=element_text(size=14),
        axis.title=element_text(size=15))

# Hides information (right pane)
mean_hiding <- svymean(~factor(gov_hiding10_V3), dat_design, na.rm = T)
confint_hiding <- confint(svymean(~factor(gov_hiding10_V3), dat_design, na.rm = T))

plot_hiding <- data.frame(x = c(0:10),
                          prop = mean_hiding[1:11],
                          lower = confint_hiding[1:11],
                          upper = confint_hiding[12:22]) %>% 
  ggplot(., aes(x = x, y = prop, ymin = lower, ymax = upper)) +
  geom_bar(stat = "identity", fill = "darkblue")+
  geom_errorbar(width=0.1, alpha=0.9, size=0.7, col = "blue2")+
  labs(x = "\n Governments are hiding information \n about the pandemic",
       y = "Percentage of respondents")+
   scale_y_continuous(labels = percent_format(accuracy = 1))+
  scale_x_continuous(breaks = 0:10, labels = c("0 \n", "1 \nExtremely unlikely", "2 \n", "3 \n", "4 \n", "5 \n", "6 \n", "7 \n", "8 \n", "9 \nExtremely likely", "10 \n"))+
   theme_minimal()+
  theme(axis.text=element_text(size=14),
        axis.title=element_text(size=15))

# Put the plots next to each other
distribution_transparency <- grid.arrange(plot_transparent, plot_hiding, ncol = 2)

ggsave(distribution_transparency, file = "./Figures/figure_A4.png", width = 10)
```

#### Figure A5: Change in the perceived probability that governments are hiding information about the pandemic between June 2020 and January 2021

```{r}
figure_A5 <- as.data.frame(svymean(~factor(gov_hiding_var23), dat_design, na.rm = T)) %>% 
  ggplot(aes(x = c(-8:6, 8, 10), y = mean, ymin = mean-1.96*SE, ymax = mean+1.96*SE)) +
  geom_bar(stat = "identity", fill = "darkblue")+
  geom_errorbar(width=0.1, alpha=0.9, size=0.7, col = "blue2")+
  labs(x = "Change in perceived probability that governments are hiding information (-10 to 10)",
       y = "Percentage of respondents")+
   scale_y_continuous(labels = percent_format(accuracy = 1))+
  scale_x_continuous(breaks = -10:10)+
   theme_minimal()+
  theme(axis.text=element_text(size=11),
        axis.title=element_text(size=12));figure_A5

ggsave(figure_A5, file = "./Figures/figure_A5.png", height = 4)
```

## Appendix B: Main results tables

#### Table B1: Factors associated with perceptions of transparency

```{r}
# Create the regression table
stargazer2(hiding_lm, transparency_lm, suppress.messages = T,
          header = F, out="./Tables/table_B1.html",
          column.labels = c("Hides information",	"Lacks transparency"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes")),
          se = starprep(hiding_lm, transparency_lm, se_type = "HC2"),
          p = starprep(hiding_lm, transparency_lm, stat = "p.value"),
          covariate.labels =c("False beliefs", "Trust scientists", "Federal PID",
                               "Provincial PID", "Ideology", 
                              "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", "Female", "Education"))
```

#### Table B2: Changes in vote intentions

```{r}
# Are false beliefs and perceptions of transparency associated with changes in vote intentions?

## Transparency perceptions and vote intentions for at the provincial level
transparency_vote_prov <- lm(vote_prov_var23 ~ gov_hiding_V2 +  trust_scien_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + region + vote_prov_V2, data = dat, weights = dat$pond_V3)

## False beliefs and vote intentions for at the provincial level
misinfo_vote_prov <- lm(vote_prov_var23 ~ misinfo +  trust_scien_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + region + vote_prov_V2, data = dat, weights = dat$pond_V3)

## Transparency perceptions and vote intentions for at the federal level
transparency_vote_fed <- lm(vote_LPC_var23 ~ gov_hiding_V2 +  trust_scien_V2 + LPC_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + region + vote_LPC_V2, data = dat, weights = dat$pond_V3)

## False beliefs and vote intentions for at the federal level
misinfo_vote_fed <- lm(vote_LPC_var23 ~ misinfo +  trust_scien_V2 + LPC_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + region + vote_LPC_V2, data = dat, weights = dat$pond_V3)

# Create the regression table
stargazer2(misinfo_vote_fed, transparency_vote_fed, 
           misinfo_vote_prov, transparency_vote_prov,
           header = F, out="./Tables/table_B2.html",
           column.labels = c("Δ Vote federal", "Δ Vote federal", 
                             "Δ Vote provincial", "Δ Vote provincial"),
           add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes", "Yes")),
           se = starprep(misinfo_vote_fed, transparency_vote_fed, 
           misinfo_vote_prov, transparency_vote_prov, se_type = "HC2"),
           p = starprep(misinfo_vote_fed, transparency_vote_fed, 
           misinfo_vote_prov, transparency_vote_prov, stat = "p.value"),
           covariate.labels =c("False beliefs", "Hides information (t-1)",
                               "Trust scientists (t-1)", "Federal PID (t-1)",
                               "Provincial PID (t-1)", "Left-right ideology (t-1)", 
                               "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", "Female", "Education", 
                               "Vote federal (t-1)", "Vote provincial (t-1)"))
```

#### Table B3: Importance of transparency perceptions

```{r}
# Estimate linear models
## Perceptions of transparency are associated with changes in trust in the federal government
transparency_trust_fed <- lm(trust_fed_var23 ~ gov_hiding_V2 + trust_scien_V2 + LPC_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + region + trust_fed_V2, data = dat, weights = dat$pond_V3)

## Perceptions of transparency are associated with changes in trust in the provincial government
transparency_trust_prov <- lm(trust_prov_var23 ~ gov_hiding_V2 + trust_scien_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + region + trust_prov_V2, data = dat, weights = dat$pond_V3)

## Perceptions of transparency are associated with changes in satisfaction with democracy
transparency_democ <- lm(satisf_democ_var23 ~ gov_hiding_V2 + trust_scien_V2 + LPC_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + satisf_democ_V2 + region + trust_fed_V2 + trust_prov_V2, data = dat, weights = dat$pond_V3)

# Create the regression table
stargazer2(transparency_trust_fed, transparency_trust_prov, transparency_democ,
          header = F, out="./Tables/table_B3.html",
          column.labels = c("Trust in the federal government", "Trust in the provincial government",	
                            "Satisfaction with democracy"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
          se = starprep(transparency_trust_fed, transparency_trust_prov, 
                        transparency_democ, se_type = "HC2"),
          p = starprep(transparency_trust_fed, transparency_trust_prov, 
                       transparency_democ, stat = "p.value"),
          covariate.labels =c("Hides information (t-1)",
                               "Trust scientists (t-1)",
                              "Federal PID (t-1)",
                               "Provincial PID (t-1)", "Ideology (t-1)", 
                              "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", 
                               "Female", "Education", 
                              "Satisfaction with democracy (t-1)",
                              "Trust federal (t-1)", 
                              "Trust provincial (t-1)"))
```

### Table B4: Factors associated with changes in perceptions of transparency

```{r}
stargazer2(hiding_var_lm, hiding_var_lm1, hiding_var_lm2, 
           header = F, out="./Tables/table_B4.html",
           column.labels = c("Δ Hides information", "Δ Hides information", "Δ Hides information"),
           add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
           se = starprep(hiding_var_lm, hiding_var_lm1, hiding_var_lm2, se_type = "HC2"),
           p = starprep(hiding_var_lm, hiding_var_lm1, hiding_var_lm2, stat = "p.value"),
           covariate.labels =c("False beliefs", "Δ Trust federal", "Δ Trust provincial",
                               "Δ Handling federal", "Δ Handling provincial",
                               "Trust scientists (t-1)", "Federal PID (t-1)",
                               "Provincial PID (t-1)", "Left-right ideology (t-1)",
                               "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", 
                               "Female", "Education", 
                               "Hides information (t-1)"))
```

#### Table B5: Factors associated with open-ended responses

```{r}
misinfo_multinom5x <- multinom(open_ended5 ~ misinfo + trust_scien_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + region, data = dat_open, maxit = 200, weights = pond_V3, Hess = T)

stargazer2(misinfo_multinom5x, 
           out="./Tables/table_B5.html",
           header = FALSE,
           column.labels = c("Handling of the pandemic", "General mistrust",
                             "Manipulation of the public"),
           add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
           covariate.labels =c("False beliefs",
                               "Trust scientists", "Ideology", 
                               "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", 
                               "Female", "Education"))
```

## Appendix C: Robustness Checks

#### Table C1: Entropy balancing for continuous treatments

```{r}
# Remove incomplete observations
dat_clean <- dat %>%
    filter(
    !(is.na(misinfo)) & 
    !(is.na(trust_scien_V3)) & 
    !(is.na(LPC_V3)) &
    !(is.na(PID_prov_V3)) &
    !(is.na(ideol_V3)) & 
    !(is.na(soc_media)) &
    !(is.na(tv)) &
    !(is.na(prints)) &  
    !(is.na(age)) &
    !(is.na(female)) &
    !(is.na(educ)) &
    !(is.na(Quebec)) &
    !(is.na(Maritimes)) &
    !(is.na(Prairies)) &
    !(is.na(BC))) 

# Estimate the entropy balancing weights
W1 <- weightit(misinfo ~ trust_scien_V3 + LPC_V3 + PID_prov_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + Quebec + Maritimes + Prairies + BC, data = dat_clean, method = "ebal", moments = 2)

# Examine whether entropy balancing successfully removed the correlations between false beliefs and the covariates
cor_before <- dat_clean %>% dplyr::select(misinfo, trust_scien_V3, LPC_V3, PID_prov_V3, ideol_V3, soc_media, tv, prints, age, female, educ, Ontario, Quebec, Maritimes, Prairies, BC) %>% 
  lapply(., cor.test, y=dat_clean$misinfo)

cor_after <- wtd.cor(cbind(W1$treat, W1$covs), weight = W1$weights)

data_balance_t <- data.frame(var = factor(c("Trust scientists", "Federal PID",
                               "Provincial PID", "Ideology", 
                              "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", "Female", "Education",
                              "Ontario", "Quebec", "Maritimes", "Prairies", 
                              "British Columbia"),
                  levels = rev(c("Trust scientists", "Federal PID",
                               "Provincial PID", "Ideology", 
                              "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", "Female", "Education",
                              "Ontario", "Quebec", "Maritimes", "Prairies", 
                              "British Columbia"))),
                        cor_before = c(cor_before$trust_scien_V3$estimate,
                                     cor_before$LPC_V3$estimate,
                                     cor_before$PID_prov_V3$estimate,
                                     cor_before$ideol_V3$estimate,
                                     cor_before$soc_media$estimate,
                                     cor_before$tv$estimate,
                                     cor_before$prints$estimate,
                                     cor_before$age$estimate,
                                     cor_before$female$estimate,
                                     cor_before$educ$estimate,
                                     cor_before$Ontario$estimate,
                                     cor_before$Quebec$estimate,
                                     cor_before$Maritimes$estimate,
                                     cor_before$Prairies$estimate,
                                     cor_before$BC$estimate),
                      p_before = c(cor_before$trust_scien_V3$p.value,
                                     cor_before$LPC_V3$p.value,
                                     cor_before$PID_prov_V3$p.value,
                                     cor_before$ideol_V3$p.value,
                                     cor_before$soc_media$p.value,
                                     cor_before$tv$p.value,
                                     cor_before$prints$p.value,
                                     cor_before$age$p.value,
                                     cor_before$female$p.value,
                                     cor_before$educ$p.value,
                                     cor_before$Ontario$p.value,
                                     cor_before$Quebec$p.value,
                                     cor_before$Maritimes$p.value,
                                     cor_before$Prairies$p.value,
                                     cor_before$BC$p.value), 
                      cor_after = cor_after$correlation[2:16],
                      p_after = cor_after$p.value[2:16]) %>% 
  gather(key = "key", value = "value", -var) %>% 
  mutate(cor = ifelse(grepl("cor", key), "Pearson's correlation", "p-value"),
         entropy = ifelse(grepl("after", key), "With entropy balancing", 
                          "Without entropy balancing"),
         entropy = factor(entropy, levels = c("Without entropy balancing", "With entropy balancing"))) %>% 
  dplyr::select(-key) %>% 
  spread(cor, value)

# Produce the balance plots
cor_plot <- ggplot(data_balance_t, aes(x = var, y = `Pearson's correlation`, 
                           col = entropy, shape = entropy))+
  geom_point(size = 2.5)+
  scale_y_continuous(breaks = c(-0.4, -0.2, 0, 0.2, 0.4))+
  scale_color_manual(values = RColorBrewer::brewer.pal(3, "Set1"))+
  coord_flip()+
  theme_minimal()+
  labs(x = "")+
  theme(axis.text=element_text(size=11),
        axis.title=element_text(size=12),
        legend.text = element_text(size=11),
        legend.title = element_blank(),
        strip.text.x = element_text(margin = margin(2, 0, 2, 0)),
        panel.border = element_rect(colour = "black", fill=NA, size=1));cor_plot

p_plot <- ggplot(data_balance_t, aes(x = var, y = `p-value`, col = entropy, shape = entropy))+
  geom_point(size = 2.5)+
  geom_hline(yintercept = c(0.05, 0.1), linetype = c(1,2))+
  scale_y_continuous(breaks = c(0, 0.1, 0.5, 1))+
  scale_color_manual(values = RColorBrewer::brewer.pal(3, "Set1"))+
  scale_x_discrete(labels = NULL)+
  coord_flip()+
  theme_minimal()+
  labs(x = "", 
       y = "p-value: Correlation test")+
  theme(axis.text=element_text(size=11),
        axis.title=element_text(size=12),
        legend.text = element_text(size=11),
        legend.title = element_blank(),
        strip.text.x = element_text(margin = margin(2, 0, 2, 0)),
        panel.border = element_rect(colour = "black", fill=NA, size=1));p_plot

# Put the plots next to each other
balance_plots <- grid_arrange_shared_legend(cor_plot, p_plot, widths = 3:2, ncol = 2,
                                            position = "top")

ggsave(balance_plots, file = "./Figures/figure_c1.png", width = 10, height = 5)

# Estimating entropy balancing with variables at t-1 (for the model measuring changes in perceptions of transparency over time)
# Remove NAs
dat_clean1 <- dat %>%
    filter(
    !(is.na(misinfo)) & 
    !(is.na(trust_scien_V2)) & 
    !(is.na(LPC_V2)) &
    !(is.na(PID_prov_V2)) &
    !(is.na(ideol_V2)) & 
    !(is.na(soc_media)) &
    !(is.na(tv)) &
    !(is.na(prints)) &  
    !(is.na(age)) &
    !(is.na(female)) &
    !(is.na(educ)) &
    !(is.na(Quebec)) &
    !(is.na(Maritimes)) &
    !(is.na(Prairies)) &
    !(is.na(BC))) 

# Entropy balancing
W2 <- weightit(misinfo ~ trust_scien_V2 + LPC_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + Quebec + Maritimes + Prairies + BC, data = dat_clean1, method = "ebal", moments = 2)

# Check balance
summary(W2)
balance_table <- cobalt::bal.tab(W2, abs = T, thresholds = 0.001)

# Estimate the models
hiding_eb <- lm(data = dat_clean, weights = W1$weights, 
   formula = formula(paste0("gov_hiding_V3 ~ misinfo +", controls)))

transparency_eb <- lm(data = dat_clean, weights = W1$weights, 
   formula = formula(paste0("transparent ~ misinfo +", controls)))

hiding_var_eb <- lm(data = dat_clean1, weights = W2$weights, 
   formula = formula(paste0("gov_hiding_var23 ~ misinfo +", controls_var)))

# Create the regression table
stargazer2(hiding_eb, transparency_eb, hiding_var_eb, 
           out="./Tables/table_C1.html",
           column.labels = c("Hides information",	"Lacks transparency", 
                             "Δ  Hides information"),
           se = starprep(hiding_eb, transparency_eb, hiding_var_eb, se_type = "HC2"),
           p = starprep(hiding_eb, transparency_eb, hiding_var_eb, stat = "p.value"),
           add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
           covariate.labels =c("False beliefs",
                               "Trust scientists", "Federal PID",
                               "Provincial PID", "Ideology",
                               "Trust scientists (t-1)", "Federal PID (t-1)",
                               "Provincial PID (t-1)", "Left-right ideology (t-1)",
                               "Social media news consumption", "TV news consumption",
                               "Print news consumption", "Age", 
                               "Female", "Education", "Hides information (t-1)"))
```

Verifying that the results are the same when including trust

```{r}
dat_clean2 <- dat %>%
    filter(
    !(is.na(misinfo)) & 
    !(is.na(trust_fed_V2)) & 
    !(is.na(trust_prov_V2)) & 
    !(is.na(trust_scien_V3)) & 
    !(is.na(LPC_V3)) &
    !(is.na(PID_prov_V3)) &
    !(is.na(ideol_V3)) & 
    !(is.na(soc_media)) &
    !(is.na(tv)) &
    !(is.na(prints)) &  
    !(is.na(age)) &
    !(is.na(female)) &
    !(is.na(educ)) &
    !(is.na(Quebec)) &
    !(is.na(Maritimes)) &
    !(is.na(Prairies)) &
    !(is.na(BC)))

dat_clean3 <- dat %>%
    filter(
    !(is.na(misinfo)) & 
    !(is.na(trust_fed_V2)) & 
    !(is.na(trust_prov_V2)) & 
    !(is.na(trust_scien_V2)) & 
    !(is.na(LPC_V2)) &
    !(is.na(PID_prov_V2)) &
    !(is.na(ideol_V2)) & 
    !(is.na(soc_media)) &
    !(is.na(tv)) &
    !(is.na(prints)) &  
    !(is.na(age)) &
    !(is.na(female)) &
    !(is.na(educ)) &
    !(is.na(Quebec)) &
    !(is.na(Maritimes)) &
    !(is.na(Prairies)) &
    !(is.na(BC))) 

W3 <- weightit(misinfo ~ trust_fed_V2 + trust_prov_V2 + trust_scien_V3 + LPC_V3 + PID_prov_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + Quebec + Maritimes + Prairies + BC, data = dat_clean2, method = "ebal", moments = 2)

W4 <- weightit(misinfo ~ trust_fed_V2 + trust_prov_V2 + trust_scien_V2 + LPC_V2 + PID_prov_V2 + ideol_V2 + soc_media + tv + prints + age + female + educ + Quebec + Maritimes + Prairies + BC, data = dat_clean3, method = "ebal", moments = 2)

cobalt::bal.tab(W3)
cobalt::bal.tab(W4)

# Estimate the models
controlsx <- "trust_fed_V2 + trust_prov_V2 + trust_scien_V3 + LPC_V3 + PID_prov_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + region"
controls_varx <- "trust_fed_V2 + trust_prov_V2 + trust_scien_V2 + LPC_V3 + PID_prov_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + gov_hiding10_V2 + region"

hiding_eb3 <- lm(data = dat_clean2, weights = W3$weights, 
   formula = formula(paste0("gov_hiding_V3 ~ misinfo +", controlsx)))

transparency_eb3 <- lm(data = dat_clean2, weights = W3$weights, 
   formula = formula(paste0("transparent ~ misinfo +", controlsx)))

hiding_var_eb3 <- lm(data = dat_clean3, weights = W4$weights, 
   formula = formula(paste0("gov_hiding_var23 ~ misinfo +", controls_varx)))
screenreg(list(hiding_eb3, transparency_eb3, hiding_var_eb3))
```

Verifying that the results are the same when using conventional entropy balancing

1. Covariate balance after entropy balancing

```{r}
# Omit NAs
dat_clean <- dat %>%
    filter(
    !(is.na(misinfo01)) & 
    !(is.na(trust_scien_V3)) & 
    !(is.na(LPC_V3)) &
    !(is.na(PID_prov_V3)) &
    !(is.na(ideol_V3)) & 
    !(is.na(soc_media)) &
    !(is.na(tv)) &
    !(is.na(prints)) &  
    !(is.na(age)) &
    !(is.na(female)) &
    !(is.na(educ)) &
    !(is.na(Quebec)) &
    !(is.na(Maritimes)) &
    !(is.na(Prairies)) &
    !(is.na(BC)) 
  )

# List the variables that you want to balance
X1 <- dat_clean[,c("trust_scien_V3", "LPC_V3", "PID_prov_V3",
                     "ideol_V3", "soc_media", "tv", "prints", "age", "female", "educ", "Quebec",
                     "Maritimes", "Prairies", "BC")]

# eBalance
eb_misinfo <- ebalance(Treatment = dat_clean$misinfo01, 
                      X = X1)

# Put the weights in a temporary data frame with unit and weight.
temp1 <- filter(dat_clean, misinfo01 == 0)
temp1$w_misinfo <- eb_misinfo$w
temp1 <- dplyr::select(temp1, record_V3, w_misinfo)

# Merge the weights and give all the treated units a weight of 1.
dat_eb <- merge(dat_clean, temp1, by = "record_V3", all.x = TRUE)
dat_eb$w_misinfo <- ifelse(is.na(dat_eb$w_misinfo),1,dat_eb$w_misinfo)

# Original treatment and control means
unwtd_means <- dat_eb %>% 
  group_by(misinfo01) %>% 
  dplyr::summarise_at(vars(trust_scien_V3, LPC_V3, PID_prov_V3,
                           ideol_V3, soc_media, tv, prints, age, female, educ,
                           Quebec, Ontario, Maritimes, Prairies, BC), mean, na.rm = TRUE) %>% 
  t() %>% 
  as.data.frame() %>% 
  slice(-1) %>% 
  rename(`Raw control` = V1, Treated = V2)

# Weighted control means
wtd_means <- dat_eb %>% 
  filter(misinfo01 == 0) %>% 
  dplyr::summarise_at(vars(trust_scien_V3, LPC_V3, PID_prov_V3,
                           ideol_V3, soc_media, tv, prints, age, female, educ,
                           Quebec, Ontario, Maritimes, Prairies, BC), ~weighted.mean(., w = w_misinfo)) %>% 
  t() %>% 
  as.data.frame() %>% 
  rename(`EB weighted` = V1)

# Create table to show balance
eb_table <- cbind(unwtd_means, wtd_means)
eb_table <- eb_table %>% 
  mutate_if(is.numeric, round, 2) %>%
  dplyr::select(Treated, `EB weighted`, `Raw control`)

rownames(eb_table) <- c("Trust scientists", "Federal ID", "Provincial ID", "Ideology","Social media news consumption", "TV news consumption", "Print news consumption", "Age", "Female", "Education", "Quebec", "Ontario", "Maritimes", "Prairies", "British Columbia")

eb_table %>% 
  kbl() %>% 
  kable_classic(html_font = "Calibri")
```

2. Results after entropy balancing

```{r}
# Estimate the models
hiding_eb_binary <- lm(data = dat_eb, weights = w_misinfo, 
   formula = formula(paste0("gov_hiding_V3 ~ misinfo01 +", controls)))

transparency_eb_binary <- lm(data = dat_eb, weights = w_misinfo, 
   formula = formula(paste0("transparent ~ misinfo01 +", controls)))

screenreg(list(hiding_eb_binary, transparency_eb_binary))
```

### Table C2: Disaggregating false beliefs

```{r}
hiding_conspi_all <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ restrict_rights_scaled + power_china_scaled + benefit_companies_scaled +", controls)))

hiding_health_all <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ co2_masks_scaled + hydroxychloroquine_scaled + false_positives_scaled +", controls)))

hiding_both_all <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ restrict_rights_scaled + power_china_scaled + benefit_companies_scaled + co2_masks_scaled + hydroxychloroquine_scaled + false_positives_scaled +", controls)))

stargazer2(hiding_conspi_all, hiding_health_all, hiding_both_all, 
          header = F, out="./Tables/table_C2.html",
          column.labels = c("Hides information",	"Hides information", 
                            "Hides information"),
          se = starprep(hiding_conspi_all, hiding_health_all, hiding_both_all, se_type = "HC2"), 
          p = starprep(hiding_conspi_all, hiding_health_all, hiding_both_all, stat = "p.value"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
          order = c("restrict_rights", "power_china", "benefit_companies", "co2_masks", "hydroxychloroquine", "false_positives", "trust_scien_V2", "LPC_V3", "PID_prov_V3", "ideol_V3", "soc_media", "tv", "prints", "age", "female", "educ"),
          covariate.labels =c("Restrict rights", "Power China", "Benefit companies",
                              "CO2 masks", "Hydroxychloroquine", "False positives",
                              "Trust scientists", "Federal PID",
                              "Provincial PID", "Ideology", 
                              "Social media news consumption", "TV news consumption", 
                              "Print news consumption", 
                              "Age", "Female", "Education"))
```

### Table C3: Disaggregating false beliefs

```{r}
hiding_conspi_var_all <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ restrict_rights_scaled + power_china_scaled + benefit_companies_scaled +", controls_var)))

hiding_health_var_all <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ co2_masks_scaled + hydroxychloroquine_scaled + false_positives_scaled +", controls_var)))

hiding_both_var_all <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ restrict_rights_scaled + power_china_scaled + benefit_companies_scaled + co2_masks_scaled + hydroxychloroquine_scaled + false_positives_scaled +", controls_var)))

stargazer2(hiding_conspi_var_all, hiding_health_var_all, hiding_both_var_all, 
          header = F, out="./Tables/table_C3.html",
          column.labels = c("Hides information",	"Hides information", 
                            "Hides information"),
          se = starprep(hiding_conspi_var_all, hiding_health_var_all, 
                        hiding_both_var_all, se_type = "HC2"), 
          p = starprep(hiding_conspi_var_all, hiding_health_var_all, 
                       hiding_both_var_all, stat = "p.value"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
          order = c("restrict_rights_scaled", "power_china_scaled", "benefit_companies_scaled", "co2_masks_scaled", "hydroxychloroquine_scaled", "false_positives_scaled", "trust_scien_V2", "LPC_V2", "PID_prov_V2", "ideol_V2", "soc_media", "tv", "prints", "age", "female", "educ", "gov_hiding10_V2"),
          covariate.labels =c("Restrict rights", "Power China", "Benefit companies",
                              "CO2 masks", "Hydroxychloroquine", "False positives",
                              "Trust scientists", "Federal PID", "Provincial PID", "Ideology", 
                              "Social media news consumption", "TV news consumption", 
                              "Print news consumption", 
                              "Age", "Female", "Education", "Hides information (t-1)"))
```

### Table C4 and C5: Exploratory factor analysis

```{r}
factor_misinfo <- psych::fa(misinfo_df, 2, n.obs = 1005, rotate="varimax")
factor_misinfo$Structure[1:6]

table_C4 <- data.frame(Statement = c("The government is exaggerating the risks\nof the coronavirus to be able to restrict\npeople’s rights and freedoms.",
                                      "The virus has been created by large\ncorporations because some of them can\ndirectly profit from it.",
                                      "The virus was created by China to increase\nits power in the world.",
                                      "The prolonged use of masks can lead to\nCO2 intoxication or oxygen deficiency.",
                                      "Hydroxychloroquine is an effective\ntreatment against COVID-19.",
                                      "Coronavirus figures are inflated because a\nsignificant number of people tested\npositive are not infected with the virus."),
                        Severity = factor_misinfo$Structure[1:6],
                        Origin = factor_misinfo$Structure[7:12],
                        Communality = factor_misinfo$communality)
rownames(table_C4) <- NULL

table_C4 %>%  
kbl() %>% 
  kable_classic(html_font = "Calibri") %>% 
  save_kable(file = "./Tables/table_C4.html", self_contained = T)
  
table_C5 <- data.frame(Factor = 1:2,
                       Eigenvalue = round(factor_misinfo$e.values[1:2], 2),
                       `Percentage of variance` = round(factor_misinfo$Vaccounted[c(2,7)],4)*100,
                       `Cumulative percentage` = round(factor_misinfo$Vaccounted[c(3,8)],4)*100)
table_C5 %>%  
kbl() %>% 
  kable_classic(html_font = "Calibri") %>% 
  save_kable(file = "./Tables/table_C5.html", self_contained = T)

dat$factor1 <- factor_misinfo$scores[,1]
dat$factor2 <- factor_misinfo$scores[,2]
```

### Table C6: Using the two factors to explain perceptions of transparency

```{r}
hiding_factor1 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ factor1 +", controls)))

hiding_factor2 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ factor2 +", controls)))

hiding_both_factors <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ factor1 + factor2 +", controls)))

stargazer2(hiding_factor1, hiding_factor2, hiding_both_factors, 
          header = F, out="./Tables/table_C6.html",
          column.labels = c("Hides information",	"Hides information", 
                            "Hides information"),
          se = starprep(hiding_factor1, hiding_factor2, hiding_both_factors, se_type = "HC2"), 
          p = starprep(hiding_factor1, hiding_factor2, hiding_both_factors, stat = "p.value"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
          covariate.labels =c("Severity factor", "Origin factor", 
                              "Trust scientists", "Federal PID",
                              "Provincial PID", "Ideology", 
                              "Social media news consumption", "TV news consumption", "Print news consumption", 
                              "Age", "Female", "Education"))
```

#### Table C7: Using the two factors to explain changes in perceptions of transparency over time

```{r}
hiding_var_factor1 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ factor1 +", controls_var)))

hiding_var_factor2 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ factor2 +", controls_var)))

hiding_var_both_factors <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ factor1 + factor2 +", controls_var)))

stargazer2(hiding_var_factor1_var, hiding_var_factor2_var, hiding_var_both_factors_var, 
          header = F, out="./Tables/table_C7.html",
          column.labels = c("Hides information",	"Hides information", 
                            "Hides information"),
          se = starprep(hiding_var_factor1, hiding_var_factor2, 
                        hiding_var_both_factors, se_type = "HC2"), 
          p = starprep(hiding_var_factor1, hiding_var_factor2, hiding_var_both_factors, stat = "p.value"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
          covariate.labels =c("Severity factor", "Origin factor", 
                              "Trust scientists (t-1)", "Federal PID (t-1)",
                              "Provincial PID (t-1)", "Ideology (t-1)", 
                              "Social media news consumption", "TV news consumption", 
                              "Print news consumption", 
                              "Age", "Female", "Education"))
```

### Table C8: Including true statements 

```{r}
hiding_susc <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_V3 ~ susceptibility +", controls)))

transparency_susc <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("transparent ~ susceptibility +", controls)))

hiding_var_susc <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding_var23 ~ susceptibility +", controls_var)))

multinom_susc <- multinom(open_ended2 ~ susceptibility + trust_scien_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ +  region, data = dat_open, maxit = 200, weights = dat_open$pond_V3, Hess = T)

stargazer2(hiding_susc, transparency_susc, hiding_var_susc, 
          header = F, out="./Tables/table_C8.html",
          column.labels = c("Hides information",	"Lacks transparency", 
                            "Δ  Hides information"),
          se = starprep(hiding_susc, transparency_susc, hiding_var_susc, se_type = "HC2"),
          p = starprep(hiding_susc, transparency_susc, hiding_var_susc, stat = "p.value"),
          add.lines = list(c("Regional fixed effects", "Yes", "Yes", "Yes")),
          covariate.labels =c("False - true beliefs",
                              "Trust scientists", "Federal PID",
                              "Provincial PID", "Ideology", 
                              "Trust scientists (t-1)", "Federal PID (t-1)",
                              "Provincial PID (t-1)", "Ideology (t-1)", 
                              "Social media news consumption", "TV news consumption", "Print news consumption",
                              "Age", "Female", "Education", 
                              "Hides information (t-1)"))

screenreg(hiding_var_susc)
```

### Figure C2: Model dependency: Density of the “False belief” estimate across 50,000 randomly sampled model specifications

```{r}
# Sensitivity Simulation
covarstosample <- c("trust_fed_V2", "trust_prov_V2", "trust_scien_V2", "LPC_V3", "CPC_V3", "NDP_V3", "PID_prov_V3", "soc_media", "tv", "prints", "radio", "age", "female", "educ", "region", "ideol_V3", "general_trust_V3", "trust_journ_V3", "RISKDOLE_V3", "satisf_life", "fear", "anger")

des1 <- svydesign(id=~1,weights=~pond_V3, data=dat)
mod1 <- svyglm(gov_hiding_V3 ~ misinfo, design = des1)


sims  <- 50000
store <- matrix(NA,sims,1)
colnames(store) <- c("OLS")
set.seed(123)
for(i in 1:sims){
 #print(i) 
        ftmp <- "gov_hiding_V3 ~ misinfo"
        bk   <- sample(covarstosample,sample(1:length(covarstosample),1),replace=F)
        ftmp <- as.formula(paste(ftmp, "+",paste(bk,collapse=" + ")))
        modols <- lm(ftmp, data = dat) 
 store[i,1] <-  c(coef(modols)["misinfo"])
}

# Estimates  
summary(store)
dx1 <- density(store[,1])

store1 <- as.data.frame(store) 
figure_C2 <-  ggplot(store1, aes(x = OLS))+
  geom_density(col = "blue3")+
    labs(x = "Estimate (false beliefs)",
         y = "Density")+
    theme_minimal();figure_C2

ggsave(figure_C2, file = "./Figures/figure_C2.png", height = 5)
```

#### Figure C3:Results of the multinomial logistic models based on the method used to deal with disagreement between coders

```{r}
misinfo_multinom <- multinom(open_ended ~ misinfo + trust_scien_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + Maritimes + Quebec + Prairies + BC, data = dat_open, maxit = 200, weights = pond_V3, Hess = T)
misinfo_multinom1 <- multinom(open_ended1 ~ misinfo + trust_scien_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + Maritimes + Quebec + Prairies + BC, data = dat_open, maxit = 200, weights = pond_V3, Hess = T)
misinfo_multinom2 <- multinom(open_ended2 ~ misinfo + trust_scien_V3 + ideol_V3 + soc_media + tv + prints + age + female + educ + Maritimes + Quebec + Prairies + BC, data = dat_open, maxit = 200, weights = pond_V3, Hess = T)

# 1. Exclude the observation
## Generate predictions
pred <- mnl_pred_ova(model = misinfo_multinom,
                      data = dat_open,
                      x = "misinfo",
                      by = 0.01666667,
                      seed = 123, 
                      nsim = 1000, 
                      probs = c(0.025, 0.975)) 

pred$plotdata <- pred$plotdata %>% 
  mutate(level = factor(open_ended, levels = c("manipulation", "mistrust", "handling", "positive"),
                        labels = c("Manipulation", "Mistrust", "Handling", "Positive")))

## Plot the results
figure_3 <- ggplot(data = pred$plotdata, aes(x = misinfo, 
                                  y = mean,
                                  ymin = lower, ymax = upper)) +
  geom_ribbon(aes(fill = level, col = level), alpha = 0.1) + 
  geom_line(aes(col = level)) +
  facet_wrap(.~ level, ncol = 2) +
  scale_color_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_fill_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  scale_x_continuous(breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  theme_bw() +
  labs(y = "Probability of falling into each category",
       x = "Average perceived truthfulness of the false statements",
       title = "Removing the observation")+
  guides(fill = "none", col = "none")+
  theme_minimal()+
  theme(axis.text=element_text(size=8),
        axis.title=element_text(size=9),
        legend.text = element_text(size=8),
         strip.text = element_text(size = 9),
        plot.title = element_text(size=9, hjust = 0.5))

# 2. Keep obs 1
## Generate predictions
pred1 <- mnl_pred_ova(model = misinfo_multinom1,
                      data = dat_open,
                      x = "misinfo",
                      by = 0.01666667,
                      seed = 123, 
                      nsim = 1000, 
                      probs = c(0.025, 0.975)) 

pred1$plotdata <- pred1$plotdata %>% 
  mutate(level = factor(open_ended1, levels = c("manipulation", "mistrust", "handling", "positive"),
                        labels = c("Manipulation", "Mistrust", "Handling", "Positive")))

## Plot the results
figure_3.1 <- ggplot(data = pred1$plotdata, aes(x = misinfo, 
                                  y = mean,
                                  ymin = lower, ymax = upper)) +
  geom_ribbon(aes(fill = level, col = level), alpha = 0.1) + 
  geom_line(aes(col = level)) +
  facet_wrap(.~ level, ncol = 2) +
  scale_color_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_fill_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  scale_x_continuous(breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  theme_bw() +
  labs(y = "Probability of falling into each category",
       x = "Average perceived truthfulness of the false statements",
       title = "Keeping coder 1")+
  guides(fill = "none", col = "none")+
  theme_minimal()+
  theme(axis.text=element_text(size=8),
        axis.title=element_text(size=9),
        legend.text = element_text(size=8),
         strip.text = element_text(size = 9),
        plot.title = element_text(size=9, hjust = 0.5))

# 3. Keep obs 2
## Generate predictions
pred2 <- mnl_pred_ova(model = misinfo_multinom2,
                      data = dat_open,
                      x = "misinfo",
                      by = 0.01666667,
                      seed = 123, 
                      nsim = 1000, 
                      probs = c(0.025, 0.975)) 

pred2$plotdata <- pred2$plotdata %>% 
  mutate(level = factor(open_ended2, levels = c("manipulation", "mistrust", "handling", "positive"),
                        labels = c("Manipulation", "Mistrust", "Handling", "Positive")))

## Plot the results
figure_3.2 <- ggplot(data = pred2$plotdata, aes(x = misinfo, 
                                  y = mean,
                                  ymin = lower, ymax = upper)) +
  geom_ribbon(aes(fill = level, col = level), alpha = 0.1) + 
  geom_line(aes(col = level)) +
  facet_wrap(.~ level, ncol = 2) +
  scale_color_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_fill_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  scale_x_continuous(breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  theme_bw() +
  labs(y = "Probability of falling into each category",
       x = "Average perceived truthfulness of the false statements",
       title = "Keeping coder 2")+
  guides(fill = "none", col = "none")+
  theme_minimal()+
  theme(axis.text=element_text(size=8),
        axis.title=element_text(size=9),
        legend.text = element_text(size=8), 
         strip.text = element_text(size = 9),
        plot.title = element_text(size=9, hjust = 0.5))

# 4. Randomly keep 1 or 2
## Plot the results
figure_3.5 <- ggplot(data = pred5$plotdata, aes(x = misinfo, 
                                  y = mean,
                                  ymin = lower, ymax = upper)) +
  geom_ribbon(aes(fill = level, col = level), alpha = 0.1) + 
  geom_line(aes(col = level)) +
  facet_wrap(.~ level, ncol = 2) +
  scale_color_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_fill_manual(values = RColorBrewer::brewer.pal(4, "Set1")[1:4])+
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  scale_x_continuous(breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  theme_bw() +
  labs(y = "Probability of falling into each category",
       x = "Average perceived truthfulness of the false statements", 
       title = "Randomly picking coder")+
  guides(fill = "none", col = "none")+
  theme_minimal()+
  theme(axis.text=element_text(size=8),
        axis.title=element_text(size=9),
        legend.text = element_text(size=8),
        strip.text = element_text(size = 9),
        plot.title = element_text(size=9, hjust = 0.5))

plots_open_ended <- (figure_3 + figure_3.5)/(figure_3.1 + figure_3.2)
ggsave(plots_open_ended, file = "./Figures/figure_C3.png", height = 6, width = 8)
```


## Other robustness checks: Alternative operationalization of the IV and DVs

#### Alternative coding of the variables

```{r}
# Perceptions of transparency
### As dummy variables
dat <- dat %>% 
  mutate(gov_hiding01_V3 = ifelse(gov_hiding10_V3 > 5, 1, 0),
         gov_hiding01_V2 = ifelse(gov_hiding10_V2 > 5, 1, 0),
### Three categories
         gov_hiding2_V3 = ifelse(gov_hiding10_V3 == 5, 0.5, gov_hiding01_V3),
         gov_hiding2_V2 = ifelse(gov_hiding10_V2 == 5, 0.5, gov_hiding01_V2))

# False beliefs
### Uncertainty (5/10) as a separate variable
dat <- dat %>% 
  mutate(restrict_rights3 = ifelse(MISINFOr1_V3 == 5, 1, 0),
         power_china3 = ifelse(MISINFOr2_V3 == 5, 1, 0),
         benefit_companies3 = ifelse(MISINFOr3_V3 == 5, 1, 0),
         co2_masks3 = ifelse(MISINFOr6_V3 == 5, 1, 0),
         hydroxychloroquine3 = ifelse(MISINFOr7_V3 == 5, 1, 0),
         false_positives3 = ifelse(MISINFOr9_V3 == 5, 1, 0),
         nunsure = restrict_rights3 + power_china3 + benefit_companies3 + co2_masks3 +
           hydroxychloroquine3 + false_positives3)
```

## Evaluate whether the results are the same when changing the operationalization of the IVs and DVs

### 1. Perceptions that governments are hiding information

```{r}
# With perceptions as a dummy and 3-point scale 
hiding_lm2 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding01_V3 ~ misinfo +", controls)))

hiding_lm3 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding2_V3 ~ misinfo +", controls)))

screenreg(list(hiding_lm, hiding_lm2, hiding_lm3))

# With the number of false beliefs endorsed as the IV
nmisinfo_hiding_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding10_V3 ~ nmisinfo +", controls)))

nmisinfo_hiding_lm2 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding01_V3 ~ nmisinfo +", controls)))

nmisinfo_hiding_lm3 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding2_V3 ~ nmisinfo +", controls)))

screenreg(list(nmisinfo_hiding_lm, nmisinfo_hiding_lm2, nmisinfo_hiding_lm3))

# Including the number of false statements that respondents are unsure about
nunsure_hiding_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding10_V3 ~ nmisinfo + nunsure +", controls)))

nunsure_hiding_lm2 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding01_V3 ~ nmisinfo + nunsure +", controls)))

nunsure_hiding_lm3 <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("gov_hiding2_V3 ~ nmisinfo + nunsure +", controls)))

screenreg(list(nunsure_hiding_lm, nunsure_hiding_lm2, nunsure_hiding_lm3))
```

### 2. Perceptions that the government lacks transparency

```{r}
# With the number of false beliefs endorsed as the IV
nmisinfo_transparency_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("transparent ~ nmisinfo +", controls)))

nunsure_transparency_lm <- lm(data = dat, weights = pond_V3, 
   formula = formula(paste0("transparent ~ nmisinfo + nunsure +", controls)))

screenreg(list(nmisinfo_transparency_lm, nunsure_transparency_lm))
```

### 3. Changes in perceptions that governments are hiding information

```{r}
# With the number of false beliefs endorsed as the IV
nmisinfo_hiding_var <- lm(data = dat, weights = pond_V3, 
  formula = formula(paste0("gov_hiding_var23 ~ nmisinfo +", controls_var)))

nunsure_hiding_var <- lm(data = dat, weights = pond_V3, 
  formula = formula(paste0("gov_hiding_var23 ~ nmisinfo + nunsure +", controls_var)))

screenreg(list(nmisinfo_hiding_var, nunsure_hiding_var))
```

### 4. Rationalization of what governments are hiding (open-ended)

```{r}
# With the number of false beliefs endorsed as the IV
nmisinfo_multinom <- multinom(data = dat_open, maxit = 200, weights = dat_open$pond_V3,
  formula = formula(paste0("open_ended ~ nmisinfo +", controls)))

nunsure_multinom <- multinom(data = dat_open, maxit = 200, weights = dat_open$pond_V3,
  formula = formula(paste0("open_ended ~ nmisinfo + nunsure +", controls)))

screenreg(list(nmisinfo_multinom, nunsure_multinom))
```